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Abstract 

In this paper we study a class of Quasi-Skipfree (QSF) processes where the transition rate 
submatrices in the skipfree direction have a column times row structure. Under homogeneity 
and irreducibility assumptions we show that the stationary distributions of these processes have 
a product form as a function of the level. For an application, we will discuss the Cox{k)/AP /l- 
queue, that can be modelled as a QSF process on a two-dimensional state space. In addition we 
study the properties of the stationary distribution and derive monotonicity of the mean number 
of the customers in the queue, their mean sojourn time and the variance as a function of k for 
fixed mean arrival rate. 
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1 Introduction 

A Quasi-Skipfree (QSF) process is a continuous time Markov process X = {Wp^O; on a two- 
dimensional state space S = {(m, i)\m £ Z, i £ {0,... where m denotes the ‘level’ of the 

state and i denotes the ‘phase’ within the level. Additionally, the jump rates are not allowed to 
cross more than one level in one direction, i.e. either the downward direction or upward direction. 
This framework is the natural extension of the embedded GI/M/1 and M/G/l-queues, and it has 
interesting structural properties in common with these processes. 

Neuts [12] investigates the embedded GI/M/l-queue as a skip-free to the right process. The matrix- 
geometric method for computing the stationary distribution of homogeneous Quasi-Birth-Death 
(QBD) processes has been applied in [12] and [11]. A homogeneous QBD process is a special case of 
homogeneous QSF processes, where the transition probabilities are not allowed to cross more than 
one level in both directions. In his book, Neuts discusses some examples of homogeneous QBD pro¬ 
cesses, e.g., the M/PH/ 1-queue and the PH/M/c- queue. He shows that the stationary distribution 
of these queues can be expressed in terms of a rate matrix. 

Recently, Katehakis and Smit [5] have introduced a new procedure to compute the invariant measure 
for a class of Markov chains. This procedure is called the successive lumping method. Further, an 
explicit solution and bounds for the steady state probabilities for the class of ergodic QSF processes 
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that possess the successive lumpability property have been derived in the paper of Katehakis, Smit, 
and Spieksnra [ID]; an analysis of the applicability requirements and numerical complexity of the 
successive lumping method is given in |9j. Ramaswami and Latouche m discuss QBD processes 
with a special structure, namely where the upward or the downward transitions rates form a row 
times column matrix. In the latter case the rate matrix can be computed explicitly. In the former 
case they show that the stationary distribution has a level product form, as is the case with the 
embedded GI/M/ 1-queue. Regarding rate matrix analysis, Ramaswami and Lucantoni, [13], exploit 
the structural properties of the G| PH) 1-queue. and extend the numerical feasibility of the matrix 
geometric approach to a wide set of problems by developing efficient schemes to compute the rate 
matrix. 

One of the main assumptions of m is that the transition rates are bounded as a function of the 
states. In the particular case of a Quasi Birth and Death (QBD) process, where the upward transition 
rates form a matrix with one non-zero row, the rate matrix cannot be computed explicitly. However, 
we will show how the stationary distribution of a higher level can be explicitly expressed in terms of 
the stationary distribution of lower levels, under an additional invertibility condition. Moreover, we 
will show that for this specific type of non-homogenous QSF processes it is possible to derive a level 
product form solution. 

We further study a phase-type inter-arrival, batch service queue, denoted PH/M Y / 1-queue, as a 
particular example of a QSF process, where the upward transition rates form a matrix with one 
non-zero row. We will consider a special case of a phase type distribution, a Coxian distribution 
(cf. [6]). Herein customers go through a maximum of k exponentially distributed phases, and after 
each phase the customer can enter the system. We will denote a Coxian distribution with k phases by 
Cox(k). Since the Coxian distributions are dense in the set of nonnegative distribution functions, we 
can approximate those by a Coxian distribution by using for example an expectation-maximization 
algorithm (see HI)- Therefore the results presented in this paper can be useful for various queues 
with different inter-arrival distributions. 

We show that the stationary distribution of the phases within a level has a product form as well, 
if the phase-type inter-arrival distribution is a Coxian distribution. In [10] the successive lumpable 
structure is specified for QSF processes, but we did not investigate the level product form that is 
derived in this paper for QSF problems with an unbounded number of levels to the left. The param¬ 
eter of the product form is the solution to a fixpoint equation, that can be numerically approximated 
efficiently. 

We conjecture that allowing for c servers will yield a product form solution based on c factors, 
similarly to the results presented by Adan et al [T[ Theorem 4.1]. Actually, we can also handle 
the Cox(k)/Ei /1 by taking the expected workload instead of the number of customers as the batch 
service distribution. 

This paper is organised as follows. In Section [2] we show that the stationary distribution of the 
levels has a product form when choosing the downward matrix D as a multiplication of a column 
vector c and a row vector r. In Section [3] we assume that there exists an exit state per level. This 
is equivalent to the existence of an entrance state per level for a revised level partition of the state 
space that keeps the QSF property intact. This means that the process is successively lumpable 
with respect to the new partition. This is used to derive an expression for the rate matrix with 
respect to the original partition in Section T3. 11 This derivation is justified by the invertibility of the 
generator of a transient Markov chain in the Appendix. Unfortunately, no explicit expression for 
the rate matrix can be derived. However, we can explicitly determine a matrix that expresses the 
stationary distribution of a lower level in terms of the stationary distribution of the higher level. 
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In Section [4] we analyze the non-homogenous Cox{k)/M 5 /1-queue, and find the parameter of the 
product form, as well as ergodicity properties. In the remainder of this section we specialise this 
queue to a queue with homogenous rates and an infinite number of phases. In Section [5] we derive 
monotonicity properties of the stationary mean number of customers in the queue and their mean 
sojourn time for fixed mean inter-arrival times. To conclude, we show that the stationary distribution 
of the -E^/M/1-queue (a special case of the Cox(k)/M Y /I where the number of phases is fixed and the 
batch size of service is 1) converges monotonically to the stationary distribution of the D/M/l-queue. 


2 The model and basic properties 


In this section we introduce the notation and derive two initial properties for the stationary distribu¬ 
tion of homogeneous QSF processes in Lemma l2.ll and Theorem 12.31 Note that by the homogeneity 
assumption the number of phases per level of the QSF process is constant, say equal to £ + 1, for 
i > 0. Without loss of generality, we may assume the QSF process to be skip-free to the left. We 
will also assume that the levels are bounded to the right, since otherwise the stationary distribution 
cannot exist, as is shown later in this section. Without loss of generality we may then assume that 
they are non-negative, i.e. m < 0. Then the infinitesimal generator Q or q-matrix (cf. [2] and [10] ) 
takes the form: 
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Elements of these 

S 


where the (£ + 1) x (£ -j 

to the s-th higher level, the same level, and the next lower level respectively, 
matrices are ufj, Wij, and dij where ufj is the (( m,i ), (m + s,j )) element from the matrix u°, Wij is 
the ((m,z), ( m,j )) element from the matrix W, and dij is the ((m,z), (m — l,j)) element from the 
matrix D , for any i,j € {0, ...,£}. The matrices U k are such that the row sum of Q is zero. For 
example these matrices can be U k ' = U l and = W + Ylili 


Throughout the paper we assume that the (/-matrix Q is irreducible, conservative, stable, and non¬ 
explosive. Additionally, we assume that the jump rates are allowed to be unbounded as a function 
of the state and that X is the minimal process. It is convenient to denote the levels by L m , so that 
L m = {(m, z) | z € {0,... , £}} and S = U m L m , where —oo < m < 0. In view of the structure of Q 
given in Eq. ()2.1I) . jumps can only take place to levels L^ for k > m — 1, given that the current level 
is Lm- 


Since the levels are mutually exclusive, they form a partition of the state space. We will introduce 
some more notation. In accordance with [lOj , for any fixed m we define the sub-level set of L m to 
be the set of states L m = U k< m Lk, while the set L m = U k>mLk is the super-level set of L m . Then 
clearly {L m _i,L m } is a partition of S for each m. 

Suppose that the QSF process is positive recurrent. Let ir denote the (unique) stationary distribution. 
We denote the m-level sub-vector of it by 7r m . The stationary distributions corresponding to L m -i 
and L m will be denoted by 7r m _i and ir m respectively. Then by a standard taboo decomposition, as 
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discussed in m, we can express 7r m _i in terms of ir m and the expected amount of time spent in the 
states of Lyii —1 • 

We denote by rnT(k,i)(l,j) the expected amount of time spent in state (l,j) without passing through 
states in the super-level set L m , given that the system starts in state (k, i), where (k, i), ( l,j) £ L m -i- 

Further, denote by „J the |L m _i| x |L m _i| matrix with elements mT(k,i)(l,j)- Write: 


0 ... 0 0 


0 ... 0 0 
0 ... 0 D 


( 2 . 2 ) 


where 0 stands for an {£ + 1) x (£ + 1) zero matrix and D_ has dimension \L m \ x \L m -i\- Then by 
the skip free property to the left 

Jim— 1 = ^m Ffrrf • ( 2 - 3 ) 

We will next express in terms of Q. The elements of the (/-matrix Q are denoted by q(k,i)(i,j): 
where (k,i), ( l,j) £ S. Denote Q(k,i) = — Q(k,i)(k,i ), the parameter of the exponential sojourn time in 
state (k,i). Further denote by m Q the taboo-generator with taboo set L m , restricted to the states 
of L m - 1 ) with elements 


_ / Q(k,i)(i,j)i (A;,i), (Z, j) £ Lm— i) 

mQ(k,i)(lj) ~ | o, otherwise, 

and m q(k,i) = ~ m Q{k,i)(k,i)- The taboo is therefore imposed at the time of the first jump out of a 
state of the QSF process on L m - 1 - By assumed irreducibility, m Q is the (/-matrix of a transient, 
non-conservative and minimal Markov process on the state space L m - 1 - Since the number of levels 
is unbounded to the left, is of infinite dimension. 

Denoting the corresponding (minimal) transition function by m P t = rnPt(ki)(lj) w here all elements 
(k,i), (l,j) £ L m - j. It follows that 


poo 

■m T (k,i)(l,j) = / mPt.,(k,i)(l,j)dt < °°- ( 2 - 4 ) 

•J 0 

We can state the following lemma. 

Lemma 2.1. The inverse matrix m Q -1 of m Q exists and m Q~ X = Further, 


mf mQ — II, Q J — b 


where \ m is the identity matrix on L m - 1 - 


Proof. By virtue of Anderson j2] (Theorem 2.2.2) the minimal transition function m P t satisfies the 
Kolmogorov forward and backward equations and it is the unique solution on S. 

Define a scalar a > 0. Then mR(a) = J 0 °° e~ at m P t dt is the associated resolvent with elements 
rrT(k,i)(l,j) (®)j where ( k,i),(l,j ) £ L m ~ i- By virtue of Anderson [2], Propositions (2.1.1) and (2.1.2), 
it holds that 

®?n^(®) — lm T mQ — lm T m/^(®) mQ' 
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Hence, by virtue of the backward equation, for (fc, i), (Z,j) G L m - 1 it holds that 


^ rnT(k,i)(l,j) (®) $(k,i)(l,j) "h mQ(k,i)(k,i) rnT(k,i)(l,j)(Q J ) 

"b ^ ^ mQ(k,i)(k',i') rnT(k',i')(l,j) (^0 5 
(k' ,i')^(k,i) 


(2.5) 


where Kronecker delta d(k,i)(i,j) 


1, ( k,i) = (l,j) 
0, otherwise. 


Observe that nf(k)(i,j)( a ) t n!T{k,i)(i,j) for a 4- 0. Using that the summation in Eq. (12.51) contains only 
non-negative terms, we can now take the limit a 0 in Eq. (12.51) and use the monotone convergence 
theorem to obtain 


\k,i)(l,j) + m Q(k,i)(k,i) rriT\(k,i)(l,j) "b E m Q(k,i)(k',i')mT(k',i')(l,j) ~ 0. 

{k',i')^{k,i) 


( 2 . 6 ) 


Eq. (12.61) is equivalent to 

^m (k,i)(l,j) T (mQ rrJ")(k,i)(l,j) 0- 

This yields that m Q = — l m . 

The relation m Q = —\ m is analogously proved, by using the forward equation for the resolvent. □ 
As a consequence of the result above we obtain: 

Em-1 = -KmUrriQ- 1 - ( 2 - 7 ) 


Define T to be the {£ + 1) x {£ + 1) subblock of (— ir Q)~ l corresponding to the states in L m -\. Note 
that without further restrictions on D it is unfortunately not possible to calculate it independently 
of the rates of states in L m - 1 - From Eq. ()2.3I) it follows that 


7T m -i = TT m B, B = DT. 

where B^j represents the expected local times spent in [m — 1 , j), before absorption into L m , given 
that the process starts in ( m,i ). 

By homogeneity, we can recursively show that 


Tim = ttqB m , for m < 0. (2-8) 

Now the existence of the stationary distribution implies that 

£ B~ m l t+1 = (/ - B)- 1 ^ < oo, 

m<0 


with I the (£ + 1) x (£ + 1) identity matrix and l^+i the i + 1-dimensional vector consisting of ones. 

The next lemma shows that it is not possible to have an unbounded state space to the right and a 
stationary distribution for a homogeneous QSF process. 

Lemma 2.2. Consider a homogeneous QSF process where the number of levels is unbounded in the 
negative direction. In order for the process to be positive recurrent, the number of levels has to be 
bounded in the positive direction. 
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Proof. Suppose that the levels are not bounded to the right. Were the stationary distribution to exist, 
the same reasoning as in the above (Eq. (12.81) 1 would yield that 7 r m = TT n B~ m+n , m < n. Analogously, 
by homogeneity it would follow that n n = ir^^a, with a an (JL + l)-dimensional positive column 
vector, independent of n (see also the proof of Theorem 12.31 belowl. Then 

^ ' 7T(m,j) 7f(n,0)®- (f -®) 1^+1 • 

m<n,j 

Taking the limit n — >• oo, on the one hand yields that Yl m <n j n (m,j) 1 - On the other hand, the 
right-hand side of the above equation is the product of two scalars, 7T(n,o) and a T (I — B)~ 1 lp + \. Since 
linvj^oo TT( n ,o) = 0, the product equals 0 as well, in the limit n —> oo. This is a contradiction. □ 


2.1 Special choice of D 


Imposing extra structure on D implies a result on the structure of the stationary distribution, that 
is described below. This extra structure covers the case when D = c ■ r, where c is a column vector 
and r is a row vector, and • denotes matrix product (of potentially non-square matrices). In other 
words, the rows of D are proportional to each other. Let the I th column of T be denoted by the 
vector U = ... ,tu) T , for i € {0 ,... ,£}. 

Then B has the following form: 


B = DT = 


c 0 r 

■ to 

c 0 r ■ 

■h • 

• • c 0 r 

■ tp 

c\r 

■ to 

cir ■ 

■ h ■ 

■ ■ c\r 

■ tp 

cpr ■ 

' to 

cpr ■ 

h ■ 

■ • cpr ■ 

■ tp 


where t = (r • to,r ■ t\,... ,r ■ tf) and 


B n = (c-t) n = c-{(i-c) n ~ 1 i) 


°i r ■ 

\i =0 


n—1 

I c-L 


n > 1 . 


Notice that matrix B has only one nonzero eigenvalue 7 = t ■ c with left eigenvector t of multiplicity 
1, provided that t ■ c is finite. Combination with Eq. (12.81) implies that the stationary distribution 
of the levels has a product form. The following theorem holds. 

Theorem 2.3. Assume that D has a column times row structure as is described above. 

If 7 < 1 and I < 00 , then the Markov process X is positive recurrent. 

If the Markov process X is positive recurrent, then 7 < 1. In either case: 

l 

71 \m,j) = 7" (m+1) Y 7r (0,i) c i h- ( 2 - 9 ) 

i =0 


Proof. Suppose that X is positive recurrent. By Eq. (12.81) . the stationary distribution necessarily has 
the form of Eq. (12.91) . Taking the summation over the levels m. in Eq. (12.91) yields a finite expression. 
This implies that E m <M < 00 for any level M. Hence, 7 < 1 necessarily. 

By separating state (0,0) from level 0, one can compute the expected sojourn time in the states of 
S \ {(0, 0)} by inverting the (negative of the) finite rate matrix restricted to this set, analogously to 
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the above derivations. Define £(o,i) to be equal to an unknown constant X( 0 ,o) times the expected 
sojourn time in state (0, i ) given a start in state (0,0) before absorption into (0,0), for all i > 0. Then 
define from Eq. (|2.9D by substituting vr( 0) j) by X( 0 ,i) in the right-hand side. Note that ®(o,o) 

is the only unknown constant. Since 7 < 1, they can be normalised to yield a (unique) probability 
distribution on S. □ 

It is still hard to compute B directly, therefore we introduce even more structure on D in the 
next section. More precisely, we assume that c has only one non-zero component. Without loss of 
generality, we may assume that cq is the only non-zero element of c. 


3 Exit states and successive lumpability 

Our main assumption concerns the presence of exit states, which are defined below. 

Definition 3.1. Let M C S. Then ( m,i ) £ M is an exit state for M, if 

i) Q(k,j)(k',j’) = 0 for all ( k,j ) £ M with (k,j) / (m,i). 

ii) 'Yh(k'Q(m,i)(k',j') > 0 . 

For the remainder of this section we make the following assumption, which implies that for all m, 
level m — 1 can only be reached directly from level m through state (m, 0). In other words, (to, 0) is 
an exit state to level m — 1 . 

Assumption 3.2. State (m, 0) is an exit state for superset L m , for all m = ..., —1,0. 


3.1 Stationary distribution 

Under assumption 13. 2\ the sub-matrix D contains only one non-zero row. This implies that matrix 
B has only one non-zero row (the first) as well. Then by virtue of Theorem 12.31 

1 ,j) (3-1) 

where = cq r ■ tj , for j = 0 ,..., £. 

Let (ij = boj/boo. We know: 7 = coto = cq r ■ to = boo < 1 and thus: 

7r (m,j)= 7~ m fijK( 0 ,o), (3-2) 

for j = 0 ,... ,£. 

As has been mentioned at the end of section [2l explicit computation of the matrix B is in general 
hard, because it can be of infinite dimension and it can have a complicated structure. In the presence 
of an exit state however, it turns out to be simpler to express n m in terms of 7 r m _i, cf. Eq. (13.31) . 

By tabooing on the set L m - 1 , there exists an |L m _i| x (£ + 1 ) matrix R m , such that 

71 m — TLm—lB m . (3.3) 

By our homogeneity assumptions, R m is in fact independent of m, and so we suppress the dependance 
on m in our further notation, whenever possible. We will next define an embedded q-rnatrix Q on L m 
and an L m - 1 x L m matrix A by removing the entrance state ( m + 1 , 0 ), such that the mean (local) 
sojourn times in the state of L m do not change, given a start in L m . We next state the following 
theorem using this shift. 
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Theorem 3.3. Consider a level homogenous process with an exit state in each level. Then: 


m — Em— 1-f^j 

where 

R ■= -AQ -\ 

where A and Q are \L m _i\ x (£ + 1) and (£+1) x. (£ + 1) matrices with elements 

oo £ , 

qij = Wij + ^ U ir ^£ ° J , » ( 3 - 4 ) 

s= 1 r=0 l^v=0 ^Ov 

oo £ , 

= u ™ : ~ k + ^2 u " ° J . ’ ^ or < m ’ 0 - ®> ^ ( 3 - 5 ) 
s=m+l—/c r=0 ^-^=0 0^ 

respectively. 

Notice that doj / ]Cr=o ^0r ^ ie probability that level set L m is entered at state (m, j), given that a 
downward transition (to set L m ) occurs starting in state (m + 1,0). 


Proof. We will explicitly compute R, containing the expected sojourn times spent in the states of 
level L m , before absorption into L m - 1 , given that the process starts in L m . 


To this end, note that since (m + 1,0) is an exit state for level L m +i, it is an entrance state for the 
set L' rn , defined below: 

L' m . = L m U {(m+ 1,0)}. 

This implies that L' m can only be reached from states in S \ L' m through (m + 1,0). The technique 
(successive lumping) developed in Katehakis and Smit [ 8 ] can be used to compute the expected 
sojourn time spent in the extended level L' m = L m U {(m + 1, 0)} before absorption into L m ~ i, given 
that the process starts at this extended level L' m . 

We use this technique to compute Q, the transition rate matrix of size (£ + 2 ) x {£ + 2 ) embedded 
on level L' m for any m. Then Q has elements: 




Wij , 

(M) 

! (^! j) ^ 4 m, 

£~1 Er= 0 <’ 

{k, i) 

G L m ,(l,j) = (m + 1,0), 

doj j 

( k,i ) 

= (m + 1,0), (7, j) € L m , 

1 

'■O. <^5 

II 

o 

a- 

o 

(M) 

= (hj) = (m + 1 , 0 ). 


(3.6) 


Note that the transitions leading to states in superlevel L m+ i\{(m + 1,0)} are mapped to the 
entrance state (m + 1,0). By virtue of Lemma 12.11 the expected sojourn time spent in level L' m 
before absorption into the sub-level set L m -1 is obtained by inverting —Q. Katehakis, Smit and 
Spieksma m, Eq. (9), show that 

T^'m = -Em-i AQ- 1 , (3.7) 

where n' m = an d the elements of A are given below with k < m: 




,m—k 
ij 5 


U. 


E oo 
S— 


m +1 


-k Ylr=0 


(U) 

(.hd) 


( m,j ) and i,j € { 0 , ...,£} 
(m + 1 , 0 ). 


(3.8) 


The validity of this exit state lumping construction of Eqs. (13. 4p and (13.51) is guaranteed by 
Lemma lA.21 in appendix A. □ 
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Restriction to QBD processes For QBD processes, where U 1 = U and U s = 0 for s > 2, the 
above derivations impy that: 

7T m = '^m— \Ri (3-9) 

where R = —UQ _1 . If U is an invertible matrix, then the result above is equivalent to 

7r m _i = - 7 TmQU- 1 . (3.10) 

In this equation we have expressed 7r m _i in terms of 7r m , similarly to Eq. (12.811 . However, QJJ~ l is 
explicitly computable, provided that U is invertible, whereas B may not be invertible. 

From Eqns. (I3.2P and (13.9p we have the following result 

(1/7)0 = W- (3-11) 

Particularly if U is an invertible matrix, then 

Pi-QW- 1 = 7 P- (3.12) 

Lemma 3.4. If l is finite, then l/y is the maximum eigenvalue of R in absolute value. 

Proof. Let R be an irreducible non-negative matrix. Then by the Perron-Frobenius theorem (for 
example in Seneta [16| ) the maximum eigenvalue is positive and real. We denote this eigenvalue by 
r. Further, there exist unique strictly positive left and right eigenvectors. Let x denote this positive 
right eigenvector. Then for alH = 0,..., t. 


rxj = y ^rjjXj. (3.13) 

j 

From Eq. (|3.11jl it is clear that I /7 is an eigenvalue of R with P, (where P > 0) as the left eigenvector 
corresponding to I/ 7 , i.e.: 

(!/ 7 ) 0 j = ^2,Pi r ij- ( 3 - 14 ) 

i 

Next, we will show that I /7 is the maximum eigenvalue of R, thus that I /7 = r. 

From Eq. (13.131) we get 

X) rXi/3i = X X PiBjXj. 

i j i 

Combining with Eq. (13.141) yields 

r(x ■ P) = — (x ■ P). 

7 

Since x ■ P > 0, it follows that r = I/ 7 . In other words, I /7 is the maximum eigenvalue of R. □ 

In the next section we will discuss a single server queueing model with Coxian arrivals and batch 
services. 
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4 Queueing application: analysis of the Cox(k)/M Y /I queue 


We consider a queueing model in which customers arrive within a maximum of k exponentially 
distributed phases. The inter-arrival distribution is a Coxian distribution of order k, see e.g. [6]. 

In a Coxian distribution of order k, phase i £ {0,..., k — 1} lasts an exponentially distributed amount 
of time with parameter A i, at the end of which either a new customer arrives with probability 1 — qi, 
or a new phase starts with probability qi, where q k -i = 0. Upon arrival of a new customer, a new 
inter-arrival distribution starts. To avoid trivialities, we assume that qi > 0 for i < k — 2. Hence the 
inter-arrival distribution is a Cox{k) distribution with parameters (Ao,..., \ k -i,qoi • • •, qk- 1 )- 

We use this naming of the phases to have a more natural state space description: here, state (m, i ) 
denotes the state of the system when there are m customers in the system and the m + 1 arriving 
customer has completed i arrival phases. We will use the random variable C ~ Cox(k ) to denote 
the inter-arrival time. 

The Cox(k)/M Y /1 queueing system can be formulated as a QSF process X{t) on the state space 
S = {Lq, Li, ...}, where L m = {{m, 0), (m, 1),..., (m, k — 1)} for all m > 0. Service occurs according 
to the distribution induced by Y ; in batches of size j, with 1 < j < b, each with probability pj. We 
will denote the probability-generating function of Y as 4>y(-). Note that ^y(x) = Yl b j=iPj x ^ ■ 

We will first analyse the general Cox(k)/M y /1-queue with finitely many phases. Then we will 
discuss some special subcases, where the exponential phases all have the same parameter, or the 
probabilities of a new inter-arrival phase are all equal upto the order. We will briefly discuss the 
extension to the case of a Coxian distribution of infinite order. 


Analysis of the Cox(k)/M Y / 1-queue The transition rate matrix Q for this model takes the 
form: 

W 0 U 0 0 0- 

D[ W U 0 0 • 

Q= D' 2 D 1 W U 0 ; ( 4-1 ) 

L>3 D 2 Di W U 


with k x k sub-matrices Dj = pi ■ pi, D' k = ^2i = i. Di, where I is the identity of size b x b and 


— Ao 

</oAo 


0 


(1 — (7o)Ao 

•• 0 

0 ' 

0 

— Ai f/i A 1 


0 


(1 - qi)\i 

0 

0 





, u = 




0 

0 

— 2 

Qk-2^k-2 


(1 — Qk- 2 )^k -2 

•• 0 

0 

0 

0 

0 

— Afc-i 


Afc_i 

•• 0 

0 


W 0 = 


and W = Wo + pi. Note that U has a c • r structure, with 

^ (1 — ^0)^0 ^ 

c = ' and r = (1,0, • • • , 0), 

(l-q k . 2 )\k -2 
V Afc-i ) 


10 











and therefore fits the analysis of Section [2j Note further that the structure of Eq. (14.11) is the 
transpose of the matrix introduced in Eq. m- We have to interchange the roles of the matrices U 
and D in the formulae of the previous sections to fit the model under consideration. Another possible 
way to look at the problem is to use the original negative numbering of the levels. Then we do not 
need the interchanging described above, but the numbering is not induced by the model anymore. 

We assume that the QSF process is ergodic. We shall prove below that for ergodicity it is necessary 
and sufficient to require that the mean number of arrivals per unit time is smaller than the mean 
number of service completions, which is specified in the relation below: 


1 


sr^k—1 
2^i =0 



Qi\ 1 


1 

EC 


b 

< *£>, = /JiEY - 

3 = 1 


(4.2) 


When ergodicity holds, we can use the results of Theorem 12.31 and the stationary distribution takes 
the form: (again, note the change of notation, since we are considering positive levels now) 


^(mj) T 


k —1 k—1 

1-1 Y = 7 m_1 Y 7r (o,i) A i( 1 “ 

2=0 2—0 


where for notational convenience, we use qk-i = 0. The factor 7 and the vector t denote the 
distribution of the phase and are implicitly given as largest positive eigenvalue and corresponding 
left eigenvector of the rate matrix R for a fixed level. Note that once 7 and t are known, the 
distribution of level 0 can be deduced. 

An alternative procedure is done by observing that the state (m, 0) is an entrance state for L m 
(and an exit state for the shifted partition L m _ 1 U {(m, 0)}). This allows us to use the results from 
Section [3] for computing the rate matrix in the reverse (downward) direction. However, we prefer not 
to use the shifted partition, because the emerging stationary distribution in both phases and levels 
will then have a notationally less amenable form, as we shall see. Therefore we will stay with the 
current level partitions. 

In consideration of the statements above, the following relation holds (cf. Eq. (13.71) 1 for levels m > 1: 

Tn — Tn-|-l2? — ^"m+1 AQ , (4.3) 

where A and Q are given in (13. 8|) and (13.61) . We refer to [10J . where a similar formula is also provided, 
but not the results regarding the product form. For this model, the matrix Q takes the form: 


Q 


—Ao qo^o 0 0 

H —(Ai + /i) qiXi ••• 0 


0 

0 


0 

jJL 0 


— (Afc-2 + At) qk-2Xk-2 

0 ~(Afc-i + At) 


(4.4) 


Note that Q is independent of the batch size distribution. The matrix A has dimension k x bk, and 
can therefore be written as: 


A = {A\ A 2 ... A/,) T , 
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where Aj is the following k x k matrix: 


Ai = 


J2 pa j ' 0 

j>i 

PiP 

j>i 

° 

j>i 


for 1 < i < b and Ab = pb • pi. 


0 

0 


PiP 


(4.5) 


Since the rate matrix R expresses the stationary solution of lower levels in terms of that of higher 
ones, we cannot recursively compute the stationary distribution in the infinite capacity case. We 
still need to find 7 and t to do so. 


Instead we will derive equations for 7 and the unknown (conditional) stationary probabilities of the 
phases. To this end we consider the expression: 

'KmQ = Tfm+lA. 

Because of the level product form solution (cf. Eq. (13.21) 1 we can rewrite this equation as: 

b b 

T^mQ — ^ ^ 'Km+jAj — TTm^^'Y Ai. (4.6) 

i=l i =1 

Using Eq. m (i-e. 7 i 7 m) j) = tj7T( m; o)) t° calculate the components of the vectors with 1 < i < k — 1 
on both sides of the above equality Eq. (14.61) yields: 


Hence 


(Ai + p)ti - qi-iXi-i ■ U- 1 = ^2 I j PjP ' U = 7) • U, 

3 = 1 

ti Qi—lXi—l 


(4.7) 


(4.8) 


U -1 Ai + p- pchil)' 

In other words, writing /% = qi-iap = U/U- 1 , we get that t t = 17;=! Ql-i a l ' to- Therefore the 


stationary distribution has the following form for m > 1 

i 

^"( m,i ) ^(1,0)7 Ql— l^Z, 

Z =1 

where a* the following function of 7 derived from (14.81) : 

Ai-i 


Oi-i — 


Xi + p - p^vi'y)' 

Note that ( 07 ,..., otk-i) only depend on (q 0 ,..., qk- 2 ) through 7 . 


(4.9) 


(4.10) 


Using the balance equation for state (m, 0) (m > 1) and Eq. (|4.10p we obtain the following expression 
for 7 as a function of (07,..., CKfc-i): 


7 = 


k —1 i 

«i(l - <?o)Ao + En am- 1(1 - 9i)Aj 

i=lZ =1 


(Ao — Ai)ai + Aq 


(4.11) 
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where we have substituted r/)y ( 7 )// = Ai + n — Ao/aq. Equations (14.1011 and (14.111) provide a system 
of k equations in the unknowns 7 and oq, z = 1,..., k — 1. As is easily checked, a* = Aj_i/Aj, 
i = 1 , — 1,7 = 1 form a solution of this system. 

Using the same balance equation, we obtain a fixpoint equation for 7 : 

7 = F( 7 ), (4.12) 


with F : R —>• R given by 


Aq( 1 ~ go) + Ei=j nU gf-i A;+/ ,-^ y(7 ) (1 - + 7M<M7) 

Aq + E 


N-i 


F( 7 ) = 


(4.13) 


Again, 7 = 1 is a fixpoint of F. The function F is tediously but easily checked to be a convex 
function of 7 on the interval [0,1 + e] for some e > 0. It is positive on this interval, and hence larger 
than the left hand side with 7 = 0 , i.e. the value of F(0). Computing the derivative at 7 = 1 yields, 
that this is strictly larger than 1 if and only if condition (14.21) holds. Hence, F has one fixpoint 
7 < 1 if and only if (14.21) holds, showing that (14.211 is necessary and sufficient for ergodicity. This 
is a standard argument used for deriving ergodicity conditions by means of a probability generating 
function approach, see e.g. m Section 4.5. 

Hence solutions 7 and {cq}A =1 can be determined (i) by determining the unique fixpoint 7 < 1 of F 
and then (ii) insert 7 into Eq. 14.101 

Next we will express the steady state probabilities of the first level in terms of the steady state 
probability in state (0,0). For sake of presentation, we omitted the derivation and present the 
results immediately: 

Ao f-7 ^ 1/n 

(414) 


and 


A° tt / t , (1 — i)4 > y{i) v'' TT 1 e ■ - 1 

*(0.0 = *( 0 . 0 )IT 11 1 + ^1 Z, 11 ai ’ for * > !’ 


s=0 


’ r ( 0 , 0 ) — 1(1 + 


7(1 - (f>Y ( 7 )) 


Ao 


k— 1 i —1 


m(1 


>+znr* 


i= 1 1=0 


-1 


(4.15) 


(4.16) 


Then by using Eqns. (14.9ft . (14.141) . (14.15ft . and (14.161) . for m> 1 we have 


Ao(l — 7)7 m 1 tt 

*( m,i ) “ ^(O.O) ,, n 11®-1«/- 


//(l-^>y(7)) ^ 


(4.17) 


We summarise all our findings above in the next theorem. 

Theorem 4.1. The Cox(k) /M* /I queue is ergodic in and only if Eq. is satisfied. If this is the 
case, the stationary distribution of the Cox(k)/M Y /I queue on levels > 1 is given by The 

factors ai, i = 1,... ,k — 1 and 7 can be calculated from Eq. \f.l2\ ) and Eq. U-10\ ). The boundary 
level can by found by Eq. H4-15\ ) and \f.l6[ ). 
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Remark 1 (Finite capacity queues). In case of a finite capacity queue of size S, the stationary 
distribution can be computed recursively, in terms of ns■ By using the fact that level S has an 
entrance state ( S , 0 ) from the right, the analysis described in the previous section can be used and 
yields that 

7r (S',i) = ~ ( 3s ,l7r (S , ,0): 


where now 


Qs 


—Ao Ao ••• 0 0 

p —(Ai + p) Ai 0 0 

P 0 ••• ~ (Afc -2 + P) 2 

p 0 • • • 0 —p 


Then using with the matrices Q m and A m associated with level m < S, one may compute 

the lower level stationary probabilities (up to a constant). Final renormalisation yields the correct 
stationary distribution. 


Remark 2 (Variable service rates). Suppose that service rates p and batch probabilities pj are equal 
to p rn and pj 1 for level m < S, j = 1 respectively, for some S < oo.Then the stationary 

distribution of Cox{k)/M 5 /I queue on the levels S + k, k > 1, can be computed in exactly the 
same manner as in the above, yielding the stationary probabilities of these levels upto a constant. 
Meanwhile, the stationary distribution for levels m < S can be computed from its by using (B 
with the matrices Q m and A m associated with level m, as described in Remark [H Again, final 
renormalisation yields the correct stationary distribution. 


4.1 A Cox(k) inter-arrival distribution with homogeneous parameters 

In this section we consider several subcases for the model described in the previous section. 


Homogeneous rates First let us assume that the rates in the exponential distribution are all 
equal: Aj = A, * = 0,..., k — 1. Then from Eq. (14.101) we get 


A 

X + p- pfirh) 


i = 1 ,..., k — 1 , 


(4.18) 


in other words, the phase factors Oj have become independent of the phase. Then Eq. (14.lip reduces 
to 

k —1 2—1 

7 = a(l — q 0 ) + ^2 ( n®) a * +1 (l “ <?*)• ( 4 -f9) 

i =1 1=0 

The function F can be simplified, but we can also directly insert Eq. (14.181) for a in Eq. (14.191) to 
obtain that the solution 7 is a fixpoint of the equation 7 = F h { 7 ) with 


F h { 7) = 


A(1 - qo ) 


k —1 2—1 


X + p- p- fiyfir) 


+ 




i =1 Z=0 


In the same manner as for the function F, we can deduce that the ergodicity condition (14.21) is 
necessary and sufficient for F h to have a fixpoint 7 < 1. 
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Homogeneous rates and probabilistic phase transitions Next we assume that additionally 
the probabilities of a new phase are all equal, that is, qi = q, i = 0,..., k — 2. This has no further 
impact on a. However, Eq. (|4.1ip reduces further to 


k—2 


1 = Y^qW +1 (l-q)+q k - 1 


- l a k . 


i=0 


The function F h becomes 
k—2 


F h h) = Y,i 


i =0 


A 


X + 


) ,+ 1 a-,)+/-/ 


A 


A + fi - /x • <Mt) 


Further, the stationary distribution of the levels rn > 1 has a product form 

A(1 - 7)7 m “ 1 (ga) < 


^( 0 , 0 ) 


and 


^(o.i) — ^(o.o )Q l 1 + 


M 1 - <M7)) 

(1 - 7 )^( 7 ) 


E 


or 


(1 - q)n(l ~ <M7)) 

where "<»'»> = (i-„*)(A + mi-<M7))) 


7(1 -</>y(7)) ^ 

(1 - g)(l - a) 


m > 1 , 


for z = 1 , • • • , A: — 1, 


1 - q k 


(4.20) 


(4.21) 


(4.22) 


Homogeneous rates and deterministic phase transitions: the Ej~/M y /1-queue Taking 
qi = q = 1, i = 0,..., k — 2, finally yields an Erlang ( k , A) inter-arrival distribution. Then Eq. (14.111) 
takes the very simple form 

7 = a k . (4.23) 


The function F h is simply given by 



_\_ 

7) 


k 


Finally, the stationary distribution has a product form (except for level 0) determined by only one 
factor 

,\(1 - 7 W+fc(m-l)) 

*<"-■> = mx-Mt)) ' m - 11 <4 ' 24) 

and 

^ = + for* = 1’-" > fc -l’ ( 425 ) 

, 1 — a 

where tt ( 0i0) = —— ■ 

Notice that this is not surprising. The process is then essentially a one-dimensional process with a 
GI/M/l-structure. This is known to have a product form stationary distribution (cf. Asmussen [3]). 

We further note that the Cox(k)/M y /1-queue with constant rates and probabilistic phase transitions 
is not of the (pure) GI/M/l-type. Indeed, the probabilities are only constant with respect to the 
phases with index smaller than k — 1 . 
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Distribution of the number of customers in the queue By taking the summation over the 
phases per level we obtain the following modified geometric distribution for the number of customers 
in Ek/M 5 / 1 -queue 


-k _ 

n m ~ 


A(1 “ 7>2 y»-‘, m > 0, 


A(1 - 7 ) 


1 - 


kn- (1 


m = 0 , 


(4.26) 


where vr^ = 7 T (rrM ), m > 0 . 


4.2 The Cox ( 00 ) /M Y /1-queue 

Allowing infinitely many phases still fits our framework. We get the natural extensions of formulas 
Eq. (14.111) and Eq. (14.131) with k = 00 . Clearly the expression for ctj as a function of 7 is not affected 
by the amount of phases. 

Restricting to the case of homogeneous rates and probabilities, so that A* = A, qt = q, i = 0,..., we 
get the following results. 

The queueing system is ergodic if and only if 

A ( < l l J2- rp y (4 - 27 ) 

v q ' 3 =1 

Then Eq. (14.111) becomes 


and so again we obtain 


7 = 


«(i - q) 

1 — qa ’ 


(4.28) 


^(7) = 


A(1 ~q) 


A(1 - q) + n(l - 4>y{i))' 


We summarise our results in the next theorem. 


Theorem 4.2. Under the ergodicity condition | J^.21 ), the stationary distribution of the Cox ( 00 ) /M^ /I 
queue on levels > 1 is given by \4-17\ ) where a and 7 can be calculated from Eq. \f.l8\) and Eq. 

The boundary level can by found by Eq. where 7T(o,o) = (1 — (/)(! — a). 


4.3 Numerical Analysis 

Assume an ergodic Cox(k)/M' 1 /I queue of finite order. In order to calculate the solutions {ai,i = 
1,..., k — 1; 7 }, one can solve (14.121) directly, e.g. by the Newton-Raphson method, and then deter¬ 
mine i = 1,..., k — 1, from (14.101) . 

As 7 is a fixpoint of the map F associated with (|4.12p . another possibility is to approximate the 
fixpoint 7 < 1 by selecting 70 and by iteratively computing 7^+1 = F( 7 „). It is simply checked that 
the fixpoint 7 = 1 is not stable, but that the fixpoint of interest, smaller than 1 , is a stable fixpoint 
. Hence we can use the following scheme to approximate the desired values oti,i = 1 ,,k — 1, 7 . 

Approximation scheme 1 
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1. Choose 7 < 1; 

2. iteratively put 7 := F( 7 ), till desired convergence; compute a t from (14.101) . i = 1, ..., k — 1. 

In the case of homogeneous rates, we have shown that the level factor 7 is the unique fixpoint smaller 
than 1 of the function F h . This leads to the following adapted scheme, that we only formulate for 
the case of homogeneous rates. 

Approximation scheme 2 for the case A* = A, i < k — 1 (where k = 00 is allowed) 

1. Choose 7 < 1; 

2. iteratively put 7 := F h ( 7 ), until desired convergence, and compute a from (14.181) . 


It is outside the scope of the paper to discuss the rate of convergence of the scheme, as well as 
a detailed stopping criterion. Table 1 below shows the non-surprising property that 7 and a are 
non-increasing in q. Further note that even for a high value of the continuation probability q , 7 is 
already approximately constant starting from around 20 phases. 


9 

k = 2 


iO 

II 


k = 20 


k = 1000 


k = 00 



7 

a 

7 

a 

7 

a 

7 

a 

7 

a 







0.1 

0.4168 

0.4414 

0.4153 

0.4411 

0.4153 

0.4411 

0.4153 

0.4411 

0.4153 

0.4411 

0.2 

0.3847 

0.4339 

0.3788 

0.4325 

0.3788 

0.4325 

0.3788 

0.4325 

0.3788 

0.4325 

0.3 

0.3538 

0.4272 

0.3406 

0.4246 

0.3406 

0.4246 

0.3406 

0.4246 

0.3406 

0.4246 

0.4 

0.3239 

0.4214 

0.3006 

0.4173 

0.3005 

0.4172 

0.3005 

0.4172 

0.3005 

0.4172 

0.5 

0.2948 

0.4163 

0.2585 

0.4105 

0.2582 

0.4104 

0.2582 

0.4104 

0.2582 

0.4104 

0.6 

0.2663 

0.4117 

0.2141 

0.4042 

0.2134 

0.4042 

0.2134 

0.4042 

0.2134 

0.4042 

0.7 

0.2385 

0.4076 

0.1673 

0.3986 

0.1658 

0.3984 

0.1658 

0.3984 

0.1658 

0.3984 

0.8 

0.2113 

0.4039 

0.1176 

0.3935 

0.1147 

0.3932 

0.1147 

0.3932 

0.1147 

0.3932 

0.9 

0.1845 

0.4006 

0.0648 

0.3890 

0.0598 

0.3886 

0.0598 

0.3886 

0.0598 

0.3886 

1 

0.1581 

0.3976 

0.0085 

0.3851 

fts 0 

0.3846 

« 0 

0.3846 

w 0 

0.3846 


Table 1. (A,/x, 7 ,pi,p 2 ,P 3 ) = (0.5,0.8,0.35,0.25,0.5,0.25). 

5 Monotonicity properties and relation with the D/M 1 /1-queue 

5.1 Monotonicity properties 

Next we will study nronotonicity properties for the homogeneous Cox(k)/M 1 /I queues, in the fol¬ 
lowing sense. If Aj = A and q- L = q for i = 0, ..., k — 2, then we denote the corresponding fc-order 
Coxian distribution by Cox(k, X,q). 

In the remainder of the paper we assume that the ergodicity condition Eq. (14.21) is satisfied. It 
is the aim to compare the stationary distribution of the number of customers in the system for 
Cox(k , A, q)/M 5 /I queues, with different inter-arrival distributions with the same mean inter-arrival 
times and associsated probabilities q, but with a different amount of phases. 
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Let A* and phase probability q and the number of phases k be given. Then, in order for the mean 
inter-arrival time to equal 1/A*, the parameter Afc of the homogeneous Coxian distribution has to be 
equal to: 


1 

A* 


1 

Afc 




1=0 


1 -q k 
Afc(l - q)' 


that is 


Afc 


A*(l — q k ) 
i -q 


Denote the corresponding factors in the stationary distribution (see previous section) of the cor¬ 
responding QSF process by ccfc and 7 fc respectively, and denote the stationary distribution of the 
number of customers in the system by Tt k . 


It does not seem possible to stochastically compare these queueing systems for a different number of 
phases: Table 2 below shows that there is a lack of monotonicity in the parameter 7 *,, especially for 
high values of the continuation parameter q. 


<7 

CM 

II 


i O 

II 


O 

II 


k = 20 


k = 50 


k = 1000 



7 

a 

7 

a 

7 

a 

7 

a 

7 

a 

7 

a 

0.1 

0.4485 

0.4734 

0.4502 

0.4764 

0.4502 

0.4764 

0.4502 

0.4764 

0.4502 

0.4764 

0.4502 

0.4764 

0.2 

0.4439 

0.4939 

0.4501 

0.5057 

0.4502 

0.5058 

0.4502 

0.5058 

0.4502 

0.5058 

0.4502 

0.5058 

0.3 

0.4371 

0.5120 

0.4494 

0.5383 

0.4502 

0.5391 

0.4502 

0.5391 

0.4502 

0.5391 

0.4502 

0.5391 

0.4 

0.4286 

0.5283 

0.4472 

0.5738 

0.4502 

0.5771 

0.4502 

0.5771 

0.4502 

0.5771 

0.4502 

0.5771 

0.5 

0.4189 

0.5429 

0.4415 

0.6111 

0.4499 

0.6206 

0.4502 

0.6209 

0.4502 

0.6209 

0.4502 

0.6209 

0.6 

0.4082 

0.5563 

0.4300 

0.6488 

0.4482 

0.6670 

0.4502 

0.6718 

0.4502 

0.6718 

0.4502 

0.6718 

0.7 

0.3968 

0.5685 

0.4105 

0.6853 

0.4413 

0.7243 

0.4499 

0.7316 

0.4502 

0.7319 

0.4502 

0.7319 

0.8 

0.3848 

0.5798 

0.3811 

0.7197 

0.4195 

0.7796 

0.4463 

0.8012 

0.4502 

0.8037 

0.4502 

0.8037 

0.9 

0.3725 

0.5902 

0.3407 

0.7514 

0.3658 

0.8307 

0.4139 

0.8738 

0.4484 

0.8905 

0.4502 

0.8912 


Table 2 . (A*, p, 'j,Pi,p 2 ,P 3 ) = (0.5,0.8,0.35,0.25,0.5,0.25). 


However, in the case of deterministic number phase transitions, an Erlang queue, there exists more 
structure regarding the change of 7 with respect to a increasing number of phases k. 

Monotonicity properties for the Ek/M Y /1-queue Let us next restrict to the case q = 1, in 
other words, the case of an Erlang inter-arrival distribution. Then the arrival rate in the /c-phase 
system is given by Afc = X*k. Write p = X* / p. 

The following results hold. 

Theorem 5.1. a) The sequence 7 fc is strictly decreasing in k. It has a limit 7 * = lim*._*.oo 7fc; which 
is the unique solution £ smaller than 1 of the equation below. 

£ = e - ( 1-< Ar(0)/p_ (5.1) 

b) The map k e-> 7 Tq = 1 — p( 1 — 7 fc)/(l — <^y(7fc)) is a strictly decreasing function if and only if 
p{y = 1 } 7 1. If P{1 — 1} — 1, i.e. the batch size equals 1 with probability 1, then 7Tq — 1 — p, 
for k = 1,2,... 
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Proof. By combination of (I4.18P and (14.231) we obtain: 


l hk = i + 


(i - <M7fc)) V 
k P ) 


(5.2) 


Define g k (x) := (l + ^ , k = 0,1,- 


Clearly this function is increasing in x. To show that g k is also increasing in k , we expand the 
expression using the binomial formula: 

k +1 


9k+ i(x) = E 


4 = 0 


k + 1 


p(Aj + 1) 


y' ( x/pY _ k k + 2 - i (x/p) k+l 
4 - ^ i! fc + 1 /c + 1 (fc + l) fc+1 


(5.3) 


On the other hand, 


sfcO) = E L 


4=0 


/c\ / x 


A>0) 


>+E 


( x/p) 1 k — 1 A; + 1 — z 


2—1 


Z! 


A; 


k 


(5.4) 


Eq. (15.31) and Eq. (15.41) yield via a term by term comparison that 

g k+ i(x) > g k (x), for all x > 0. (5.5) 

In other words, g k+ i(l - <Mt)) > g k { 1 - (pyh)), with g k+1 (l - 4>y{ 1)) = 3fc(l - <Ml)) = 1 and 
g k+ i(l — 4>y(0)) = (1 + p _1 ) k+l > (1 + p^) k ■ Recall that 7 j is the unique hxpoint of the equation 
1/7 = 5i(l-^y ( 7 ))) with 7* e (0,1), * = k, k + 1. For x < 7 *, 1/x > ^(l-0y(x)) and for x € ( 7 *, 1), 
necessarily 1/x < gi( 1 — fy(x)), i = k, k + 1. Hence 7^+1 < 7 &. 

Since 7 ^ are non-increasing, and bounded below, the sequence has a limit, 7 * say, with 7 * < 1 . This 
limit solves (15.11) by the standard limiting argument that lim*,_>.oo(l + x/k) k = e x . The function 
7 1 —y Q-P-^Ypi)) i s a convex function, with hxpoint s 7 = 1 and 7 * < 1 , derivative larger than 1 at 
7 = 1, and positive value at 7 = 0. The result then follows in a standard manner, thus completing 
the proof of a). 

Part b) follows from the fact that 

1 — 7 1 

1 -^( 7 ) = eC«£7W 

which is strictly increasing in 7 < 1 if and only if P{E = 1} < 1. 

Clearly, with increasing k, the variance of the inter-arrival time decreases. Hence, the average server 
utilisation strictly improves with decreasing variance, in the case of batch sizes Y, with P {Y = 1} < 1. 
Whereas, if the batch size equals 1 with probability 1, the average server utilisation is equal to p and 
thus constant. 

We can say a little more. We define L k = m ■ Tr k r to be the mean number of customers in the 
system under the stationary distribution. Further, we denote by W k the expected sojourn time and 
by 14 the variance. The following comparison result holds. 


(5.6) 

□ 
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Theorem 5.2. The following are true: 


Lk +1 < L k W k +1 < Wk, and, V k+1 < V k , for k= 1,2,.... (5.7) 

Proof. The expectations with respect to 7 ffc_|_i and n k are equal to L k+ \ and L k respectively. 

It is easy to check that 

Lk = i- T~t —7’ for k = l, 2 ,... . (5.8) 

1 - <m7fc) 

Since 7 k are non-increasing, then 


1 - <M7fc) < 1 - <M7fc+i)- 


It follows that L k+ 1 < L k . Application of Little’s formula yields W k+ \ < W k . 

From (14.261) it is clear that the number of customers in the E k /M/ 1-queue is (almost) geometrically 
distributed. This implies that V k = p 2 'y k /(l — 'y k ) 2 - It is easy to check that this yields V k+ \ < V k . 

Next, we need to check the variance for a general batch size distribution: 


V k 


P(1 -lk) 2 

1 - <M lk) 


E “hr 1 

m>l 



Applying geometric series sum yields 

E^hr 1 = 7» E - W 2 + E “T*" 1 

m>l m>l m>l 

27fc 1 

(l- 7 fc ) 3 (1 — 7fc ) 2 ' 


(5.9) 


(5.10) 


By combination of (15.81) . (15.91) . and (15.101) 

= p{ l-7fc ) 2 f ^lk 1 \ _ P 2 

1 - <M7fc) U 1 - 7fc ) 3 (1 ~ lk) 2 J (1 - <M7fc )) 2 

P f 1 + 7 fc _ P \ 

1 - (j) Y (j k ) V 1 - 7fc 1 - 0Y(7fc )) ' 


We define /i( 7 ) := 


1 + 7 _ P 

1-7 1-<M7) 


Since 7 *, > 7 * and 7 ^ | 7 * by virtue of Theorem 15.II lab it is sufficient to show that h'{ 7 ) 
for 7 € [ 7 *, 1 ). 


h'{ 7) = 


(1-7 ) 2 


2-p 




Applying Eq. (15.61) yields 


h'( 7) = 


(1 -7) S 


2-P 




2 6 
i=i 


dfeW >Q 

<17 - 


(5.12) 


20 
























Since j 7 J 1 < Y^iJa 7* f° r 7 ^ [ 7 *, 1), it follows that 


h'{ 7 ) > 


1 

(1^7? 



For 7 <G [ 7 *, 1), 

_ P _< _ P _ 

E$=i^EI=o 7' " E$=i^EU(r) 1 ' 

So, to show that h'( 7 > 0) for 7 £ [ 7 *, 1 ), it is sufficient to show that 

2 >_ P _=_ P _ 

' E‘.i Pi EfcoV)‘ E‘„ ftW ’ 

or 

2(l-^( 7 *))/p>l-7*. (5.13) 

By virtue of Theorem 15.II (a), 7 * > 1 — (1 — <?V( 7 *))/p (since e~ x > 1 — x, for x > 0 ). 

This implies that 

(1 - <M7*))/P > 1 -7*- (5-14) 


So that (15.131) follows. 

Then 

T^i(7) ,i(7t) - I - M- IM ) hbt+l) - 

In other words, we have proved that 14 > Vk+i- □ 

Interestingly enough, for general batch size distribution, the stationary distribution does not stochas¬ 
tically decrease with increasing k. This follows immediately from the fact that Y m >i * s strictly 
increasing in k, whereas 7 *, is strictly decreasing. However, if the batch size is indentically equal to 
1, then the stationary distribution has a stochastically monotonic behaviour as a function of k. 

Corollary 5.3. Suppose that P{Y = 1} = 1. The following is true for k = 1,2,...: 

St 

TT/c+1 — 

or equivalently: Y m >M 7 ^ +1 < E m>M ^m> f or al1 M = 0,1,.... 

In this particular case we can also prove that a*, is strictly increasing in k. 

Lemma 5.4. Suppose that P {Y = 1} = 1. The sequence of parameters a.k is strictly increasing in 
k. 


Proof. From Eq. (14.181) we have a & 


A k 

A k + n- /ry k ' 


We define 

fk{x) = -^x k+1 - (1 + y~)x + 1. 

Afc Afc 


(5.15) 
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Taking the first derivative of /*. yields 




Then at = 
Thus 


/ i+Afc \ 1//fc 

a* k = ( (k+i) ) P°i n f where the polynomial fk{x) has a unique minimum. 


otk £ 


/i+hfc\ 1//fc / \ 

For any k > 1, we have a* k = ( J = (^1 - ^-(1 - p)J 


l + t 


Oi* k 


We define 


= ( 1 - FTT (1 “ P) 


1/k 


1/k 


We will show that g(k ) is stricly increasing in k. 

g'(k) = g(k ) (-p log (l - p^-(l - p)) - 


1 ~P 


k(k + !)(! + kp) 


Since g(k) = a* k > 0, it is sufficient to show that 

1 / k . 

-F k>g T^TT ll ■'’ , 


1 ~ P 


k(k + !)(! + kp) 


> 0. 


We know that: 


This implies 


1 


> 


1 


k + 1 (A: + !)(! + kp) 


, for k > 1. 


1-P 


1-P 


6 (fc+l)(l+A:p) g fc+1 ]> j 1 — 


fc A 1/fc 

(1-P) 


k + 1 




1-P 


fe(/c + !)(! + fcp) 


(5.16) 


(5.17) 


(5.18) 


(5.19) 


(5.20) 


(5.21) 


By combining (15.201) and (15.211) it is clear that g(k) is stricly increasing in k. This means that 
a* k < a* k+l and so a* k+l "f 1, for k —> oo. 

Next, for any positive k and x £ (0,1), the functions fk and fk+\ have 2 intersection points: at 0 
and 1. This follows from the relation: 


fk+i(x) - f k (x) = 


1 


k(k + 1 )p 


x{\ — x) 2 ^2 ' lxl 1 • 


2—1 


This also implies that fk+\{x) > fk(x) for all x £ (0,1). Hence < ak+ 1 < 1. 

In other words we can say that ak is strictly increasing in k and the proof is complete. 


(5.22) 


□ 
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5.2 Comparison of batch service queues with Erlang arrivals versus a determin¬ 
istic inter-arrival time 


Let us assume for the moment that the batch size is identically equal to 1, i.e. P{T = 1}. Consider 
the D/M/ 1-queue with mean inter-arrival time 1/A* and mean service time 1/p, where A* = A k/k. 
In this case it is quite well-known that the stationary distribution of the E^/M/ 1-queue, with the 
same mean inter-arrival and mean service time distribution, converges setwise and weakly to the 
stationary distribution of the D/M/ 1-queue, as k —> oo. We will generalise these results when there 
are batch services. By virtue of Asmussen [3] and Bhat |5] the stationary distribution i to of the 
D/M/ 1-queue is given by 


f (1 — a)pa m 1 , m> 0 
l 1 - P, m = 0, 


where a is the unique root smaller than 1 of 


a = e -^~ a)lp , 


(5.23) 


when p < 1. By virtue of Theorem 15.11 this means that a = 7 * = lim^oo and in particular 
a < 7 ^.. It follows directly that: 

7f D 5 vr fc+1 A 7f fc , A: = 1,2,... (5.24) 


Define Lo, Wo, and Vo as the mean number of customers, the expected sojourn time, and the 
variance of customers of the D/M/ 1-queue respectively, under the stationary distribution. We now 
derive the result below. 

Corollary 5.5. The following are true: 


i) Eq. \5.24\) holds for all k = 1,..., and 

ii) Lk f Lo, WkfWo, and 14 i Vo, as k —> 00 , monotonically. 

Intuitively it seems clear that a similar result holds as well for the E^/M Y /l- and D/M Y /1-queues, 
with the same mean inter-arrival times and the same batch service distributions. So far, we have 
not been able to find any result on (setwise and weak) convergence of the stationary distribution 
of the Ek/M Y /1-queue to the stationary distribution of the D/M^ /1-queue, although it should be 
completely similar to convergence results in the case of batch size equal to 1 . 


General batch size Suppose now again that P{y = 1} < 1. 

Theorem 5.6. For the stationary distribution tt d of the D/M ^ /1-queue the following holds: 


i) n!/ n —>• 7 r^, for m = 0,... and 


ii) 


7t£! 


p(l - a ) 2 1 

l-^y(cr) 

Pi 1 ~ cr ) 

1 - 4>y{a)’ 


m > 0 , 


1 
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m = 0 


(5.25) 







where a is the unique root smaller than 1 of 


b 

a = e _ ( 1_ ^ >y ( <T ))/ p , </>y(cr) = ^^PjCP. (5.26) 

3 = 1 

Proof. Since the proof is quite standard, we only provide a sketch of the proof. 

First, notice that the embedded processes on arrival instants are of the GI/M/l-type considered in 
[K pp. 82-86]. As has been derived there, it follows for the corresponding stationary distribution 
if// of the D/M 1 /1-queue that 

*A,m = I 1 - 

with a the unique root smaller than 1 of Eq. (15.261) . 

Secondly, by using semi-regeneration, cf. Asmussen (3, Ch. VII.5 and p. 283], we then obtain 
formula Eq. (15.251) for the stationary distribution of the non-embedded D/M 1 /1-queue. The de¬ 
sired convergence properties then follow from the explicit formulae for the respective stationary 
distributions. □ 

It directly follows from Theorem l5.6l il and Theorem 15.21 that the result of Corollary 15.51 holds in the 
batch service case as well. 

Corollary 5.7. The monotonicity result in Corollarv \5.5\ (ii) holds in the case of batch service. 

As a consequence, the mean number of customers, expected sojourn time and variance are minimized 
by deterministic inter-arrival times. This is a well-known result for non-batch systems (cf. Asmussen 
[3] pp. 336-339]). Results on other performance measures are given in this section as well. 
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A Appendix 

Lemma A.l. Let X be a possibly non-conservative, non-explosive, transient, stable Markov process 
on state space S with q-matrix Q. Then Q^ 1 exists and Q~ 1 = T where the ( i,j)-th element of T is 
defined as follows: r M - = Jq° P t ,u,j)dt < oo, where Pt,(i,j) are the elements of the transition function 

Pt- 

Proof. A generalization of Lemma 12.11 and the proof is analogous. □ 

Lemma A.2. Let X be a possibly non-conservative, non-explosive, transient, stable Markov process 
on state space S with q-matrix Q. Let s £ S be given. Consider the transition rate matrix Q on 
A\{s} where the elements are given by: 


Qij 


Qij T qis 


Qsj_ 
ds 1 
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where q s = -q ss . 


Let q - 1 , q - 1 denote the ( i,j)-th elements of matrices Q 1 and Q 1 ; respectively. Then the following 
are true for i,j ^ s: 

• Qij — Qij ) 

• a -1 = V M 1 

q „ ^rj ■ 


Proof. By Lemma IA.1I the inverse matrices Q -1 and Q _1 exist, and the entries are equal to the 
expected sojourn time spent in each state of S and S s = S\{s} respectively. It is convenient to use 
a representation based on the jump chain on S with transition matrix P J . The entries of P J are 
given by: 


J 

Pi,j = —: 

J q% 


(A.l) 


where % = -q iti , for i, j € S. 


Clearly P J is a sub-stochastic matrix, since Q is transient. We denote its ?r-th iterate by P J ' n . The 
0-th iterate is the identity. It follows from Anderson |2], Proposition (4.1.1) (see also Eq. 112.411 of 
this paper) that: 


T i,j — 


n>0 


J,n 1 


Qj 


(A.2) 


The jump transition probability matrix P J associated with Q on S s , is given by 


Pi . j Pi,j Pi .sPs.j ■ 


(A.3) 


This is precisely the transition matrix of the jump chain associated with Q, embedded on S s . By 
using Eq. (IA.2I) and the fact that the (i,j)-th element of —Q _1 represents the expected amount of 
time spent in state j, given a start in state i before absorption outside S, we directly obtain that 
— Q~ l is the expected amount of time spent in S s before absorption outside this set. 

This proves the first statement. For the second statement, we invoke (IA.2I) . Then for j ^ s 


‘s,3 


^ZPs.r T r,J ~ Q qs Tr ’ s ' 

r^s r^s 


The result follows. 


□ 
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